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Table of contents 


@ ¢ HISTORY ; Detailed Current Edit Histor 

(3) ATIONS — ; Declarative Part of odule 

(4) 1 § MTHSSQRTR2 = Standard Single Precision Floating SQRT 
(5) 1 MTHSSQRT_R2 = JSB SQRT routine 


984 AX/VMS Macro V04-00 P 
g-SEP 984 ¢: 98 8 ERTHRIL SRC JUVXSQRTR2.MAR; 1 
-TITLE UVXS$SQRTR2 s penn Point Square Root routine 
IDENT /1-016/ : File: MTHSQRTR2.MAR EDIT: JCW1016 


(AA RARRAAAARRAALALAAEALALALASASE SESE RAE SESE SESE ERE RSE R EE SE EERE REESE EET ESS SS 


COPYRIGHT (c) 1978, 1980, 1982, 4 BY 
DIGITAL EQUIPMENT CORPORATION, MAYNARD. MASSACHUSETTS. 
ALL RIGHTS RESERVED. 


THIS ope: I$ uit UNDER A LICENSE AND MAY BE USED AND Sd a 
ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE 
INCLUSION OF THE ABOVE Copy cHt NOTICE. THIS SOFTWARE OR ANY OTHER 
COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY 
Paanisteheree’ NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY 


THE ne eat IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE 


AN H BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT 
CORPORATION. 


DIGITAL ASSUMES NO RESPONS . 


IBILITY 
SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. 


:* 
:* 
-* 
-f 
-f* 
-® 
-* 
-® 
-® 
»2 
-* 
-® 
-* 
-* 
-® 
-* 
-* 
-® 
-* 
-* 
-* 
-* 


FOR THE USE OR _ RELIABILITY OF ITS 


LARRRRBAAAAAAAAALEALARASLALELALASSSAS ELSES EERE EE RARER ERR RARER ESR RRR RES ARES SS 


FACILE: MATH LIBRARY 
“ABSTRACT: 
MTHS$SSQRT_R2 is a special routine which is the same as MTH$SQRT except 


a faster non-standard JSB call is used with the argument in RO and no 
registers are saved. 


SONA WNE WN HO ODNAMNE WIN OD OONAUES WN | O ODNAU EW OWOONOULS Wh 


heehee ay once nb Sicilia ims? cle ies Bina -Rislog Bete Ge Ge Ge Ge Ge Ge Ge Ge Ge Ge Ge Se Ge Ge Ge Ge Ge Ge Ge Ge Ge 


UES REAR EES PWWWWIIIIIinononononononononony 3 2 Ss Se 


SOOCOOOOOSOOOOSOOOSOOSOOSOOSOSOOOO OOOO OSOOOOOCOOOOOOCOO 


SOOO OCOOOOOOOOOOCOoOoOoOSooSoSo 


VERSION: 01 
HISTORY: 
AUTHOR: 
Peter Yuo, 15-Oct-76: Version 01 


MODIFIED BY: 


Ol-1 Peter Yuo, 22-May-77 
Q1-2 Peter Yuo, 31-May- 
01-16 Jeffrey C. Wiener 14-Jul-83 
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 petailed current Edit History 6-SEP-1984 1 9% 98: 33 EMTHRTL. SREJOVXSORTR2.MAR: 1 " 3) 


-SBTTL HISTORY ; Detailed Current Edit History 


3; Flo 
HISTO 


ALGORITHMIC DIFFERENCES FROM FP=11/C ROUTINE: none 
Edit History for Version 01 of MTH$SQRT 


te Code saving after code revi 
1-2 ROTL sh ft in iSn afte into INghest bit. Use ASHL instead. 
oo ineteyess on after ADJU Py anes? .§ changed. ipte ADDW to prevent 
or flow i :16> = PFFF an 16> = 
i i nish Bik Phenat ing TOrjuneci 977 
hb gery hm A to »  eorresnesi 
eee Cc an e 
Changed otrar banal ten mechanism. Put error result in RO before 
aes’ ng MTHSS$SIGNAL e order to ie user modify error result. 


gi7$ Return -0.0 on negati ve arg. TNH 20-Dec- 
1-7 Edit in Rich Lary’s code bums. JSB routine is now _R2. JMT 19-Jan-78 
01-9 Move .ENTRY symbol to module header. TNH 14- ug: -73" 
1-010 - Put version number in Bg ye format: three ight edit 
numbers. Also ypdet nay ht notice. JBS 16-NOV-78 
1-011 - change, ATH SQUROONEG to TaTHSR Qu OONE 6; eclBs S 07-DEC-78 
ae $4 - o the PSECT directive. JBS -DEC- 
1-3) - Declare oo s. SBL 17-May-1979 
1-014 = Move MTH$SQRT_R2 to separate 7 ea var on 
pF - Change external references to te s 6-0ct-8 


Changed the reference to a D_ tloating. Be wa to a G_floatin 
instruction. Also made other minor changes in the code to accomodate 
this change. This was done to conform to ae | policy that G_floating, 
ys b float ing. should be used to backup F_floating point foutines. 


SOCGOCOOCOOOOCOOOCOOOCOSOOOOOOOOSOOSSOSoSO *7e 
SOSCDOCOOCOSOSOOOOOCSOSOOSOOOOOCOOOSOOSOSOOO «<- 
SOOoOoooooooooooooooooooooQoooooooo 
SOoooooooooooooooooooooooooooooo 


NO OONAUE WN O ODNAUE WIN  O OONOAUE wry 
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6 
XSSQRTR2 Floating Point Square Root ro = 16-SEP-1984 :08: AX/VMS Ma v04-00 p 
wie bet taating Pong Square Root routine sur GrsEpctone 7:58.38 EATANTE caelovKsoRsR2.man;1 2% 
¢ ~SBTTL DECLARATIONS 3; Declarative Part of Module 
; INCLUDE FILES: 
0; 
4! : EXTERNAL SYMBOLS: 
38 ' .DSABL GBL 
0 94 «EXTRN MTHSK_SQUROONEG 
6 32 sEXTRN MTHSSSIGNAL 
00 3 : 
4 98 ; EQUATED SYMBOLS: 
00 99 ; 
4 100 ; MACROS: none 
00 +4 3 
444 19 3 PSECT DECLARATIONS: 
00000000 104 -PSECT _MTHSCODE PIC, SAR,LONG, EXE ,NOWRT 
0000 105 3 program section for math routines 
0000 196 ; 
0000 107 ; OWN STORAGE: none 
0000 108 ; 
44 109 ; CONSTANTS: 
000 110 
0000 117 3 
43 i : Constants A and B chosen for k = odd 
13CD5FD4 0000 114° LF_ODD_A_E63_ = “X13CD5FD4 
3€4A2018 Bons 12 LF-ODD-B_EM63 = *x3C4A2018 
38 11? : Constants A and B chosen for k = even 
F61A4015 0000 119 ° LF_EVEN_A = “XF61A4015 
4B8231FD7 0000 120 LF-EVEN_B_EM64 = *X4B231FD7 
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X$SQRTR2 ; Floating Point Square Root routine 16-SEP-1984 02:08: AX/VMS Macro V04-00 p 6 
016 ATHSSORTRE "= Standard Single Precision G-SEP=1984 11:29:28 UMTHRIL.SRETOVXSGRTRZ.MAR:1 2° (a) 


-SBTTL MTHSSQRTR2 = Standard Single Precision Floating SQRT 


+0 | 
; FUNCTIONAL DESCRIPTION: 
: SQRT = single precision floating point function 
SQRT(X) is computed using the following approximation technique: 
If X <® 0, error. Let X = IXi. 
Let X = 2e*K * F where F is the fractional part. 


If K = even, X = 2**(2P) « F 
SQRT(X) = 2eeP * SORT(F), 1/2 =< F< 1 


If K = odd, X = 2"8(2P41) ® F = 208(2P42) ® (F/2) 
SQRT(X) = 2e8(P41) ® SQRT(F/2), 1/4 =< F/2 < 1/2. 


Let F* = AeF + B, 


A = 0.453730314(octal), 
B = 0.327226214(octal), for K = even. 
= Ae(F/2) + B, 
A= 8 Boor ehtoctal;. 
4 B = 0.230170444(octal), for K = odd. 
an 
K' =P, for K = even 
=P +1 for K = odd. 


Let YO = 2**k' * F* as a staight line ging wthin the 
given interval using coefficients A and B which minimize the 
absolute error at the midpoint and endpoint. 
Starting with YO, two Newton-Raphson iterations are performed. 
YCn+1] = (1/2) * ( YEn) + X/Y¥EnJ]) 
The relative error is < 10**-8, 
} 


; CALLING SEQUENCE: 


0 sart.wf.v = MTHS$SQRTR2(x.rf.r) 

33 ; INPUT PARAMETERS: 
00000004 0000 LONG = 4 ; define Longuorg multiplier 
00000004 0 x = 1 * LONG ; Contents of x is the argument 


IMPLICIT INPUTS: none 
OUTPUT PARAMETERS: 

VALUE: floating square root of the argument 
IMPLICIT OUTPUTS: none 


ee ed ad td ad = dd dS 8 8 8 8 a 2 do 
NNN NAA AAAAAAAAUNIUUINUIUIUIUIUTE BE BEE EE EE WIN IInonononununony 
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Standard Single Precision 6=SEP-1984 11:29: MTHRTL.SRCJUVXSQRTR2.MAR; 1 (4) 
9 


COMPLETION CODES: none 

SIDE EFFECTS: 

Signals: MTHS_SQUROONEG if X < 7° with reserved sgerens in RO (copied to 
the signal mechanism vector CHFSL_MCH_RO/R1 oy et AR SIGNAL) 


Associated message is: ‘SQUARE ROOT OF NEGAT VALUE". Result is reserved 
operand -0.0 unless a user supplied (or any) error handler changes CHFSL_MCH_RO/R1 


yogssen tne 
-016 


=z 
— 
zm 
Ann 
“ro 
oo 
nao 
wv 
'o 


NOTE: This procedure disables floating point underflow, enables integer 
overflow, causes no floating overflow or other arithmetic traps, an 
preserves enables across the call. 


kk kk 8 
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-SBTTL MTHSSQRT_R2 = JSB SORT routine 
; JSB SQRT = used by the standard, and directly. 


CALLING SEQUENCE: 
pove anything in RO:R2 


eee R 
JSB mTH$sQRT_R2 


ng Point $ 
Ro = isB” 


: input in RO 


Note: This routine is written to avoid any integer overflows, floating overflows, 
floating underflows or divide by 0 conditions, whether enabled or not. 


GOODOOCOOOCOOOnO 


3; second iteration, do in double precision to get truncated( rather than 
3; rounded) result. 


| 
} 7 
; 
4 
; return with result in RO 
bp 10 ; REGISTERS USED: 
0 11 RO - Floating argument then result 
00 \§ Ril = X saved for use during iteration 
44 1? R2 - scratc 
00 15 MTH$SQRT_R2:: ; JSB routi for SQRT 
i. a 099 16 Move RO, R1 3; test +) > af X and save it in R1. 
sp 15 903 1 BLEQ ZERO_NEG ; branch te ZERO_NEG if X =< 0 
005 219: x>0 
88538 fos 
005 ¢ 3 CLRL ~(SP) z make room for 2nd half of 
005 3 :; double length operand 
52 50 3¢€ 0005 4 MOVZWL RO, R2 ; isolate low 16 bits Hi exp,>fract) in R 
5 94 0008 5 CLRB 3 R2 now has Aly and left ? exp bits 
=~ 3 AA QOOA § BICW R2, RO ; clear sign and left 7 exp bits 
5 95 000D TSTB R 3: check low bit of exp 
10 18 OOOF 8 BGEQ EVEN 3; and branch if 1 
13CD5SFD4 BF 44 0011 9 MULF #LF_ODD_A_E63, RO 3; add 64 (half of bias) to (exponent-2) 
B38 0 3 and start approximation calc 
3C4A2018 8F 40 001 1 ADDF #LF_ODD_B_EM63, RO ; RO = (first approx) * 2**-64 
ii oot $ BRB ADJOST 3 go adjust 
8 1 4 EVEN: 
0 2000 4 AO 0021 5 ADOW #*x2000, RO ; exp is 0 - make it 64 (2**-64) for legalit 
F61A401 F 44 88 6 § MULF #LF_EVEN_A, RO 
4B231FD7 8F 40 ° ? on nest ADDF #LF_EVEN_B_EM64, RO ; RO = (first approx) * 2**-64 
i. a wee $0 4 33 ” ROTL #31, R2, R2 3 djviée Re Coxpep ios) by 2, 
: giving (ex 
50 52 a0 a} ADDW R2, RO : Taser exer? in first approx and 
bee ¢g 3; re-bias it. 
038 rt: :; first iteration - single precision is sufficient 
52 51 50 47 46 3 DIVF3 RO, R1, R2 ; R2 = x/Y0 
0 2 + F 4 ADDF _R ; RO = YO + Hyg 
50 080 BF A $ re SUBW #*x80, RO : RO = Y1 = (1/2)(¥O + X/Y0) 
4 4 3 no overflow possible 
; 
+ $88 


——_—_—_—_—} 
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XSSQRTR2 ; Floating Point Square Root ro rine 16-SEP-1984 AX/VMS v04 P 
oie ATHSSORT AZ = JSB SORT routine. ~3p=198¢ 17:99:38 EMTMRTE. exesovKSaRTR2.maR;1 2% 5) 
47 S 333 CLRL Re : lower part (xX) = 0 
47 5 333 DIvd RO, R1 3 div ge Y1 into X with low-order 
47 § ; 32 bits of Y1 garbage. This doesn't 
of : the df accuracy, since Y1 innacurate 
3 anyw 
47 § : PUSHL R11 3 convert x and place on stack 
47 60 ; CLRL R1 3; clear low pert of v1 
47 61; DIVD3 RO, (SP)+, R1 3; divide Y1 int 
rf $6 3 DF R1, RO ; RO = Y1 + higher. part(Xx/Y1) 
a er 
re $5 8 EHR Sau See BiYeg ADL Spe, cs 
50 ef $3 bee 98 CVTGF RO, R 
50 0080 BF A osc % SUBW #*x80, RO : RO = SQRT(X) = (1/2) (Y1 * X/Y¥1) 
0 bees 4 SQRTX: RSB 3: return, RO = result 
ee Hpi rc 
$b 75 2ERO_NEG: 
13 006 6 BEQL SQRTX 3: return with RO = result = 0 
6—E DD 0064 7 PUSHL (SP) : return PC from JSB routine 
7E OO'8F 9A 0066 78 MOVZBL #MTHSK_SQUROONEG, -(SP) ; condition value 
50 01 OF 78 006A 79 SHL #15, #T, RO : = result = reserved operand -0.0 
43 80 3 0 99 es to ° signet mechanism vector 
06E 81 3 RCH $L_MCH_RO/R1) so error handler 
006E 3 3 con modi iy the result. 
00000000°GF 02 FB 0075 o CALLS #2, G*MTHSSSIGNAL : signal over ond vee reek user" s PC 
: ependent o 
05 007 85 RSB : joeee = - RO restored from CHFSL_MCH_RO/R1 
Bee seo 
076 87 - END 
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Symbol table 


a 
1 


SEP=1986 14:09:28 UMIHATL. sae 


$ 
—J 
= 

‘ 
vr 
o 
Ww 
nnununn 


geereree 


L 
MTHSS$S1GNAL 


K SQUROONEG *eeeeeee 


ae 8 
ae RT_R2 8888000 r } 
1 


dovooe2 


SQRTX 
ZERO_NEG 


PSECT name Allocation PSECT No. Attributes 
Ss 00000000 ( 0.) 00 ( 0.) NOPIC USR CON ABS LCL 
MTHSCODE 00000076 ¢ 118.) 01 ¢ 1.) PIC USR CON’ REL LCL 
ww oonm ence rete nw emo mma 
H Performance indicators ! 
Phase Page faults CPU Time Elapsed Time 
Initialization 31 00:00:00.08 00:00:01.39 
Command pig 116 00:00:00.64 Bo BR 8° £0 
Pass 1 79 00:00:00.79 8208 383-30 
TS homespun of HEARS O 88:88 
Symbol table output 00: 0:00.02 00:00:00:02 
Psect synopsis output 00:00:00.0 00:00:00.1 
Cross-reference output 9 sb 64 es Be So 
Assembler run totals 29 0:00:02. 0:00:12.16 


The working eee Limit was 900 pages. 

3796 bytes (8 pages) of virtual memory were used to buffer the intermediate co 

There tere 10 pages of symbol table space allocated to hold 14 non-local and 8° local symbols. 
source Lines were read in Pass 1, producing 8 object records in Pass 2 

1 page of virtual memory was used to "define 1 macro. 


Serene rem nam neem mares somes > 


Macro Library name 
_$255$DUA28: CSYSLIBISTARLET.MLB;2 0 
0 GETS were required to define 0 macros. 


There were no errors, warnings or information messages. 


0 Vv04-00 
JUVXSORTR2. MAR; 1 


a => NORD NOWRT NOVEC BYTE 


RD NOWRT NOVEC LONG 


med, 


Lee 


>>>> 


>> 


tO = POM OO 
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MACRO/ENABLE=SUPPRESSION/DISABLE=(GLOBAL , TRACEBACK) /LIS=LIS$:UVXSQRTR2/O0BJ=0BJ$:UVXSQRTR2 MSRCS$:MTHJACKET/UPDATE=(ENHS:MTHJACKET) +MS 
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